View the html here: https://raw.githack.com/jess-devine/BIO4000W_GIS/refs/heads/main/README.html
Is Southern Afrotemperate Forest Expanding in the Cape Peninsula?
Investigating whether pioneer forest species indicates forest expansion. Using GIS and species occurrence data to assess forest dynamics.
Species that are Southern Afrotempraet Forest pioneers: Virgilia oroboides (Cape Keurboom), Virgilia divaricata (Garden Route Keurboom), Rapanea melanophloeos (Cape Beech), and Kiggelaria africana (Wild Peach).
Observation data from iNaturalist and City of Cape Town vegetation layers will be used.
## Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
##
## Attaching package: 'gridExtra'
##
##
## The following object is masked from 'package:dplyr':
##
## combine
##
##
##
## Attaching package: 'cowplot'
##
##
## The following object is masked from 'package:patchwork':
##
## align_plots
##
##
## The following object is masked from 'package:lubridate':
##
## stamp
##
##
## Loading required package: viridisLite
species_list <- c("Virgilia oroboides", "Virgilia divaricata",
"Rapanea melanophloeos", "Kiggelaria africana")
species_plots <- list()
species_data <- list()
# get veg data
veg <- st_read("data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous.shp")
## Reading layer `Vegetation_Indigenous' from data source
## `/home/jess/GIT/BIO4000W_GIS/data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1325 features and 5 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -63972.95 ymin: -3803535 xmax: 430.8125 ymax: -3705149
## Projected CRS: WGS_1984_Transverse_Mercator
st_write(veg, "data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_duplicate.shp", append = FALSE)
## Deleting layer `Vegetation_Indigenous_duplicate' using driver `ESRI Shapefile'
## Writing layer `Vegetation_Indigenous_duplicate' to data source
## `data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_duplicate.shp' using driver `ESRI Shapefile'
## Writing 1325 features with 5 fields and geometry type Polygon.
vegr <- st_read("data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_Remnants.shp")
## Reading layer `Vegetation_Indigenous_Remnants' from data source
## `/home/jess/GIT/BIO4000W_GIS/data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_Remnants.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 3428 features and 7 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -63951.23 ymin: -3803532 xmax: 420.7595 ymax: -3705506
## Projected CRS: WGS_1984_Transverse_Mercator
st_write(vegr, "data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_duplicate.shp", append = FALSE)
## Deleting layer `Vegetation_Indigenous_duplicate' using driver `ESRI Shapefile'
## Writing layer `Vegetation_Indigenous_duplicate' to data source
## `data/cape_peninsula/cape_peninsula/veg/Vegetation_Indigenous_duplicate.shp' using driver `ESRI Shapefile'
## Writing 3428 features with 7 fields and geometry type Polygon.
# Crop vegr
ext <- c(-66642.18, -3809853.29, -44412.18, -3750723.29)
names(ext) <- c("xmin", "ymin", "xmax", "ymax")
vegr <- st_crop(vegr, ext)
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
for (species_name in species_list) {
# Retrieve iNaturalist data
inat_data <- get_inat_obs(taxon_name = species_name,
bounds = c(-35, 18, -33.5, 18.5),
maxresults = 1000)
# Filter iNaturalist data
inat_data <- inat_data %>%
filter(positional_accuracy < 46 &
latitude < 0 &
!is.na(latitude) &
captive_cultivated == "false" &
quality_grade == "research")
# Convert to spatial object
inat_data <- st_as_sf(inat_data, coords = c("longitude", "latitude"), crs = 4326)
inat_data <- st_transform(inat_data, st_crs(vegr))
inat_data <- st_intersection(inat_data, vegr) # Remove points in urban areas
# Ensure National_ is a character (fix color mapping issue)
inat_data$National_ <- as.character(inat_data$National_)
# Store in the list
species_data[[species_name]] <- inat_data
# buffer points not in Southern Afrotemprate Forest
nsaf <- inat_data %>%
filter(National_ != "Southern Afrotemperate Forest") %>%
st_buffer(dist = 250)
#Intersect new polygons with veg remnants and filter for those that overlap Southern Afrotemperate Forest only
nsaf <- st_intersection(nsaf, vegr) %>% filter(National_.1 == "Southern Afrotemperate Forest")
# Get count of non-buffered points per vegetation type
no_buffer_summary <- inat_data %>%
st_drop_geometry() %>%
group_by(National_) %>%
summarise(count_no_buffer = n())
# Get count of buffered points per vegetation type
buffer_change_summary <- nsaf %>%
st_drop_geometry() %>%
group_by(National_) %>%
summarise(count_buffer_change = n()) # Renamed to avoid confusion
# Compute total buffered points
total_buffered <- sum(buffer_change_summary$count_buffer_change, na.rm = TRUE)
# Merge summaries
combined_summary <- full_join(no_buffer_summary, buffer_change_summary, by = "National_")
combined_summary <- combined_summary %>%
mutate(count_no_buffer = replace_na(count_no_buffer, 0),
count_buffer_change = replace_na(count_buffer_change, 0))
# Adjust count_buffer correctly
combined_summary <- combined_summary %>%
mutate(count_buffer = ifelse(National_ == "Southern Afrotemperate Forest",
count_no_buffer + total_buffered, # SAF gets the initial count plus total buffer count
pmax(0,count_no_buffer - count_buffer_change))) # Others get adjusted values
# Rename columns for clarity
combined_summary <- combined_summary %>%
rename(
"Vegetation Type" = National_,
"Initial Count" = count_no_buffer,
"Buffered Change" = count_buffer_change,
"Final Count" = count_buffer
)
# Create the Plot
species_plot <- ggplot() +
annotation_map_tile(type = "osm", progress = "none") +
geom_sf(data = inat_data, aes(color = National_), shape = 16) + # Uniform shape
labs(title = paste("Observations of", species_name)) +
theme_minimal()
# Create the Table
table_plot <- tableGrob(combined_summary)
# Print Separately
print(species_plot) # Print the figure
# Move to a new page & draw the table
grid.newpage()
grid.draw(table_plot)
}
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
All of the tree species, except Rapanea melanphloeos, had more observations more than 250m from Southern Afrotemprate Forest than within this vegetation type. However, the extent of Southern Afrotemprate Forest is small on the Cape Peninsula.
The majority of species seem to be expanding into Peninsula Sandstone Fynbos more than other vegetation types. However, Kigillaria africana seems to be expanding more into Peninsula Granite Fynbos.
# Combine all species data into a single dataframe
combined_data <- bind_rows(species_data, .id = "Species")
# Define a color palette for vegetation type
veg_colors <- viridis::viridis(length(unique(combined_data$National_)))
# Define unique shapes for each species
species_shapes <- c("Virgilia oroboides" = 16, # Filled circle
"Virgilia divaricata" = 17, # Filled triangle
"Rapanea melanophloeos" = 18, # Filled diamond
"Kiggelaria africana" = 15) # Filled square
# Generate the final plot
ggplot() +
annotation_map_tile(type = "osm", progress = "none") +
geom_sf(data = combined_data,
aes(color = National_, shape = Species),
size = 3) + # Adjust size for better visibility
scale_color_manual(values = veg_colors) + # Assign vegetation colors
scale_shape_manual(values = species_shapes) + # Assign species shapes
labs(title = "All Four Species Observations",
color = "National Vegetation Type",
shape = "Species") +
theme_minimal() +
theme(legend.position = "right")
# Ensure data is in WGS84 (longitude-latitude)
combined_data <- combined_data %>% st_transform(4326)
vegr <- st_transform(vegr, 4326)
# Define a color palette for speceis
library(RColorBrewer)
species_palette <- colorFactor(brewer.pal(min(9, length(unique(combined_data$Species))), "Set1"), combined_data$Species)
# Define a color palette for vegetation type
veg_palette <- colorFactor(viridis(length(unique(vegr$National_))), vegr$National_)
# Create iNaturalist links for popups
combined_data <- combined_data %>%
mutate(iNat_link = paste0("<a href='https://www.inaturalist.org/observations/", id,
"' target='_blank'>View on iNaturalist</a>"))
# Leaflet Map
leaflet(combined_data) %>%
addTiles() %>%
addPolygons(
data = vegr,
fillColor = ~veg_palette(National_), # Use colorFactor function
color = ~veg_palette(National_),
weight = 1,
fillOpacity = 0.5,
popup = ~paste("<b>Vegetation Type:</b>", National_),
group = "Vegetation Remnants"
) %>%
addCircleMarkers(
data = combined_data,
lng = ~st_coordinates(geometry)[,1],
lat = ~st_coordinates(geometry)[,2],
radius = 5,
color = ~species_palette(Species), # Use colorFactor function
fillOpacity = 0.9,
popup = ~paste("<b>Species:</b>", Species, "<br>",
"<b>Vegetation Type:</b>", National_, "<br>",
iNat_link),
group = "Observations"
) %>%
addLegend("bottomright",
pal = species_palette,
values = combined_data$Species,
title = "Species") %>%
addLegend("topright",
pal = veg_palette,
values = vegr$National_,
title = "Vegetation Type") %>%
addLayersControl(
overlayGroups = c("Vegetation Remnants", "Observations"),
options = layersControlOptions(collapsed = FALSE)
)
From these maps it apears that Southern Afrotemprate forest may be expanding in the Northern half of the Cape Peninsula. Historical pine plantations that are now nature reserves in this region may be hotspots for forest expansion, likely because fire was excluded for a long period of time while the plantations were in operation. Additionally, Silvermine dam may be facilitating Southern Afrotemprate Forest expansion.
However, some of the observations are in parks, botanical gardens and plantations. This likely accounts for the majority of points in Cape Flats Sand Fynbos. Aditionally, it is likely natural for there to be Southern Afrotemprate Forest to exit around small streams and rivers but the vegetation map scal is too course to map these.